%%
%This file extracts data from harmonics analysis block from PSCAD.
% Parameters needs to be set as follows in Harmonic Impedance Component
% Frequency Increment Type: Linear Scale
% Impedance output type: Phase Impedance
% Frequency output unit: Hz
% Impedance output units: Magnitude and Phase (Deg)


% Created By Farzeen Faiz Adi Rajah
% For any clarifications email:  Farzeenfaiz15@gmail.com

%% Data Extraction
clear
a=fopen('Filename.out')             %Replace Filename with Output file name from PSCAD
ref=[0,0,0,0,0,0,0,0,0,0,0,0,0];
text1=fgetl(a)
 while feof(a)~=1
     text1=fgetl(a);
     text2=fgetl(a);
     text3=fgetl(a);
     if isempty(text1)||isempty(text2)
         break
     end
     c1=textscan(text1,'%f64');
     c2=textscan(text2,'%f64');
     c3=textscan(text3,'%f64');
     num1=c1{1};
     num1=num1';
     num2=c2{1};
     num2=num2';
     num3=c3{1};
     num3=num3';
     num=[num1,num2,num3];
     ref=[ref;num];
 end
 ref(1:2,:)=[];
 plot_data=ref;
 fclose('all')
 %% Plotting the data
 Freq_increment=1;     % Enter Freq.Increment value from PSCAD
 %These numbers below maps the coloums to corresponding impedances
 Mag_Z11=2;      Angle_Z11=3;
 Mag_Z22=6;      Angle_Z22=7;
 Mag_Z33=12;     Angle_Z33=13;
 Mag_Z12=4;      Angle_Z12=5;
 Mag_Z13=8;      Angle_Z13=9;
 Mag_Z23=10;     Angle_Z23=11;
 Freq_axis=Freq_increment*plot_data(:,1);
 
 %% Magnitude plots
 figure('Name',"Magnitude plots")
 subplot(3,2,1)
 plot(Freq_axis,plot_data(:, Mag_Z11))
 title('Z11');hold on;grid minor;
 
 subplot(3,2,2)
 plot(Freq_axis,plot_data(:, Mag_Z12))
 title('Z21');hold on;grid minor;
 
 subplot(3,2,3)
 plot(Freq_axis,plot_data(:, Mag_Z22))
 title('Z22');hold on;grid minor;

 subplot(3,2,4)
 plot(Freq_axis,plot_data(:, Mag_Z13))
 title('Z31');hold on;grid minor;

 subplot(3,2,6)
 plot(Freq_axis,plot_data(:, Mag_Z23))
 title('Z32');hold on;grid minor;
 
 subplot(3,2,5)
 plot(Freq_axis,plot_data(:, Mag_Z33))
 title('Z33');hold on;grid minor;
 
 %% Angle Plots
 figure('Name',"Angle plots")
 subplot(3,2,1)
 plot(Freq_axis,plot_data(:, Angle_Z11))
 title('Z11');hold on;grid minor;
 
 subplot(3,2,2)
 plot(Freq_axis,plot_data(:, Angle_Z12))
 title('Z21');hold on;grid minor;
 
 subplot(3,2,3)
 plot(Freq_axis,plot_data(:, Angle_Z22))
 title('Z22');hold on;grid minor;

 subplot(3,2,4)
 plot(Freq_axis,plot_data(:, Angle_Z13))
 title('Z31');hold on;grid minor;

 subplot(3,2,6)
 plot(Freq_axis,plot_data(:, Angle_Z23))
 title('Z32');hold on;grid minor;
 
 subplot(3,2,5)
 plot(Freq_axis,plot_data(:, Angle_Z33))
 title('Z33');hold on;grid minor;